\tema{M'etodos iterativos para sistemas lineales}

\subsection{Repaso}

\begin{itemize}
	\item Sea $A \in \R{n}{m}$, $\lambda \in \Re$ es autovalor si $\exists x \neq 0$ 
				tal que $Ax = \lambda x$. $\lambda$ es autovalor y $x$ es un autovector.
				Para calcularlos, calculo los ceros del polinomio caracter'istico: 
				$det(A - \lambda I)$
				
	\item \texttt{Radio espectral de a:} $e(A) = max |\lambda|$, donde $\lambda$ 
				es autovalor de A. Adem'as $e(A) \leq ||A||$ (norma inducida).
				
	\item Sea $A \in \R{n}{n}, e(A) < 1 \Rightarrow (I - A)$ es no singular y 
				$(I - A)^{-1} = \sumatoria{k = 0}{\infty} A^{k}$
				
	\item Decimos que A es convergente si $\lim_{k \rightarrow \infty } (A^{k})_{ij} = 0$
	
	\item A es convergente $\Leftrightarrow$ $e(A) < 1$
\end{itemize}

\intro

Sea el sistema $Ax = b$. Queremos encontrar una sucesi'on tal que 
$\{\x{k}\} \rightarrow_{k \rightarrow \infty}  x^{*}$ \fl

\subsection{M'etodos}

En esta secci'on se usar'a la notacion $\x{ij}$ donde $i$ es la posici'on y $j$ el 
n'umero de iteraci'on.

\subsubsection{Jacobi}

\Formula{
\x{i k+1} = \frac { b_{i} - \sumatoria{j = 1}{i - 1} a_{ij} \x{j k} - \sumatoria{j = i+1}{n} a_{ij} \x{jk} }
									{a_{ii}}
}

\FormaMatricial{
Se descompone A de la siguiente manera:
A = matriz con los elementos de la diagonal (D) - matriz triangular inferior con los elementos por debajo 
de la diagonal (con el signo cambiado) - matriz triangular superior con los elementos por encima
de la diagonal (con el signo cambiado).

De esta manera, se obtiene la forma matricial:
\[
Ax = b
\]

\[
(D - L - U) x = b
\]

\[
Dx - (L+U)x = b
\]

\[
Dx = b + (L + U)x
\]

\[
x_{k+1} = D^{-1}b + D^{-1} (L + U)x_{k}
\]
}

\subsubsection{Gauss-Seidel}

Un an'alisis sobre el m'etodo de Jacobi sugiere una mejora al algoritmo.
Si queremos calcular $x_{i}^{k}$ utilizamos las componentes de $x^{k-1}$. 
Como para $i > 1$, ya se calcularon $x_{1}^{k}, \hdots, x_{i-1}^{k}$, y 
probablemente sean mejores aproximaciones de las soluciones reales 
$x_{1}, \hdots, x_{i-1}$ que $x_{1}^{k-1}, \hdots, x_{i-1}^{k-1}$, parece
m'as razonable calcular $x_{i}^{k}$, por medio de los valores calculados
m'as recientemente. Esto es el coraz'on de Gauss-Seidel y se muestra en la 
siguiente f'ormula:

\Formula{
\x{i k+1} = \frac { b_{i} - \sumatoria{j = 1}{i - 1} a_{ij} \x{j k+1} - \sumatoria{j = i+1}{n} a_{ij} \x{jk} }
									{a_{ii}}
}

\FormaMatricial{
Se descompone A de la siguiente manera:
A = matriz con los elementos de la diagonal (D) - matriz triangular inferior con los elementos por debajo 
de la diagonal (con el signo cambiado) - matriz triangular superior con los elementos por encima
de la diagonal (con el signo cambiado).

De esta manera, se obtiene la forma matricial:
\[
Ax = b
\]

\[
(D - L)x - Ux = b
\]

\[
(D - L)x = b + Ux
\]

\[
x_{k+1} = (D - L)^{-1}b + (D - L)^{-1} Ux_{k}
\]
}

\typ

Las formas matriciales de Jacobi y Gauss-Seidel responden a $x_{k+1} = c +Tx_{k}$ \fl

\prop{
La sucesi'on $\{ x_{k} \} $ definida por: $x_{k+1} = c +Tx_{k}$ converge a $x^{*}$
soluci'on de $x = c +Tx$ $\Leftrightarrow$ $e(T) < 1$.
} \fl

\prop{
Si A es diagonal dominante $\Rightarrow$ Jacobi converge.
} \fl

\prop{
Si A es diagonal dominante $\Rightarrow$ Gauss-Seidel converge.
} \fl

\prop{
Sea $||T|| < 1$ norma inducida, entonces
\begin{itemize}
	\item $\x{k} = Tx_{k-1} + c$ converge
	\item $||x - \x{k}|| \leq ||T^{k}|| ||\x{0} - x ||$
	\item $||x - \x{k}|| \leq \frac{||T^{k}||}{1-||T||} ||\x{1} - x ||$
\end{itemize}
}

\subsubsection*{Error en m'etodos iterativos}

Desde le punto de vista intuitivo paree razonable que, si $\tilde{x}$ es una 
aproximaci'on a la soluci'on x de $Ax = b$ y si el vectorial residual $r = b - A\tilde{x}$
tiene la propiedad de que $||r|| = ||b - A\tilde{x}||$ es peque'no, entonces $||x - \tilde{x}||$
tambi'en ser'a peque'no. A menudo este es el caso, pero algunos sistemas, que en 
la pr'actica se presentan con frecuencia, no poseen esta propiedad.

\prop{
Supongamos que $\tilde{x}$ es una aproximaci'on a la soluci'on de $Ax = b$, que
A es una matriz no singular y que r es el vectorial residual de $\tilde{x}$. 
Entonces, para toda norma natural,
\[
||x - \tilde{x}|| \leq ||r|| ||A^{-1}||
\]
}

\subsection{Gradientes Conjugados}

\subsubsection*{Introducci'on}

El m'etodo del gradiente conjugado de Hestenes y Stiefel fue desarrollado 
originalmente como un m'etodo directo dise'nado para reoslver un sistema lineal 
$n \times n$ definido positivo. Como m'etodo directo, por lo general es inferior 
a la eliminaci'on gaussiana con pivoteo, pues ambos m'etodos requieren n pasos para 
determinar una soluci'on y los pasos del m'etodo del gradiente conjungado tienen 
un costo mayor en c'alculos que los de la eliminaci'on gaussiana.

Sin embargo, el m'etodo del gradiente conjugado es muy 'util como m'etodo iterativo 
de aproximaci'on para resolver sistemas esparcidos de gran tama'no con entradas 
no nulas que aparecen con patrones predecibles. Cuando la matriz est'a precondicionada 
para que los c'alculos sean m'as eficaces, se obtienen buenos resultados aproximadamente
en $\sqrt{n}$ pasos. Empleado de esta manera, este m'etodo es preferible sobre la
eliminaci'on gaussiana y los m'etodos iterativos ya analizados.

\subsubsection*{Desarrollo}
Consideremos el sistema $Ax = b$ y busquemos la soluci'on pensandola como el
m'inimo de una funci'on $Q(x): \Re^{n} \rightarrow \Re$.

Sea $Q(x) = x^{t}Ax - 2x^{t}b$. Para hayar el m'inimo de esta funci'on debemos 
analizar la derivada primera (y estudiar los puntos cr'iticos, es decir, donde
es igual a 0); y la derivada segunda (para ver que efectivamente es un m'inimo).

Para esto tenemos que:

\[
\nabla Q(x) = 2 (Ax - b) = 0
\]

Como se puede observar, los puntos cr'iticos de Q se dan cuando se resuelve el 
sistema para el que se buscaba soluci'on.

Para ver que es m'inimo, derivamos de vuelta y vemos que:

\[
\nabla^{2} Q(x) = 2A
\]

\todo{C'omo sigue esto????}

Ahora que vimos que esta forma de resolver un sistema de ecuaciones, veamos 
c'omo encontramos el m'inimo efectivamente.

Para esto definimos la siguiente sucesi'on:

\[
\x{k+1} = \x{k} + \alpha_{k} + d_{k}
\]

El objetivo ahora es obtener una forma de encontrar $\alpha_{k}$ y $d_{k}$.

Se puede ver facilmente que:

\[
\alpha_{k+1} = \frac{\alpha_{k} (b-A\x{k}) } {\alpha_{k}^{t} A d_{k}}
\]

Veamos c'omo se obtienen las direcciones conjugadas:

\subsubsection*{Direcciones conjugadas}

\defi{
Sea A sim'etrica definida positiva, decimos que $u_{1}, \hdots, u_{n}$ son 
A-conjugadas si:
\[
u_{i}^{t}Au_{j} = 0 \forall i \neq j
\]

Si adem'as, $u_{i}^{t}Au_{i} = 1 $, se llaman A-ortonormales.
} \fl

\prop{
Sean $d_{1}, \hdots, d_{n}$ direcciones A-ortonormales.
Sea $\x{k} = \x{k-1} + \alpha_{k-1} + d_{k-1}$, entonces $x_{n}$ es la soluci'on 
de $Ax = b$ para cualquier $\x{0}$.
} \fl

\prop{
Sea $r_{k} = b -Ax_{k}$, entonces $r_{k}^{t}d_{i} = 0 \forall i = 0, \hdots, k-1$
} \fl

Luego, para generar las direcciones tenemos las siguientes f'ormulas:

\[
d_{k} = - r_{k} + \beta_{k} d_{k-1}
\]

donde:

\[
r_{k} = b - A\x{k}
\]
\[
\beta_{k} = \frac{r_{k}^{t} A d_{k-1}} {d_{k-1}^{t}Ad_{k-1}}
\]

utilizando $d_{0} = - r_{0}$. \fl

\prop{
\begin{itemize}
	\item $<r_{0}, r_{1}, \hdots, r_{k}> = <r_{0}, Ar_{0}, \hdots, A^{k}r_{0}>$ (<>: subespacion generador por)
	\item $<d_{0}, r_{1}, \hdots, d_{k}> = <r_{0}, Ar_{0}, \hdots, A^{k}r_{0}>$
	\item $d_{k}^{t} A d_{i} = 0 \forall i = 0, \hdots, k-1$
\end{itemize}
}